Methods and systems for determining heave and heave rate of vessels

ABSTRACT

A method for determining heave and heave rate for a vessel is described. The vessel includes an inertial navigation system (INS), and sensors for the INS being located at a point A of the vessel, the vessel having a zero heave reference point B, and the described method provides heave and heave rate at a point C of the vessel. The method includes determining reference coordinates for points A, B, and C of the vessel, generating a velocity signal representative of a velocity at a point on the vessel, generating a heave rate signal based upon the velocity signal, and generating a heave based upon the heave rate signal.

BACKGROUND OF THE INVENTION

[0001] This invention relates generally to navigation of vessels, andmore specifically, to determination of and compensation for heave andheave rates for vessels.

[0002] Heave is the vertical distance of motion relative to sea-level.and a vertical motion of a vessel riding waves at sea, relative to sealevel, is generally referred to as heave rate. For a vessel riding thewaves, the vertical motion is simply an oscillation about sea-level,having a zero long-term average of earth-relative displacement andvelocity. Heave and heave rate values are utilized in controllingvarious operations, such as, to control a tether to a submerged diver,submarine, and for other underwater operations, as well as to define amotion of a sonar transmitter, and to define deck motion at a point forlanding operations of a helicopter.

[0003] Some known inertial systems are used to measure heave and heaverate. However, a vertical channel of a pure inertial system typicallyhas divergent errors. A global positioning satellite (GPS) system can beused to compensate for such divergence, but the GPS system itself mayintroduce errors into the heave and heave rate determinations.

BRIEF SUMMARY OF THE INVENTION

[0004] In one aspect, a method for determining heave and heave rate fora vessel is provided. In an exemplary embodiment, the vessel is equippedwith an inertial navigation system (INS), and sensors for the INS arelocated at a selected point of the vessel, e.g. point A. The vesselfurther has a zero heave reference point B. The method is utilized todetermine heave and heave rate at a point C of the vessel, where point Cis a different point than points A and B. More particularly, the methodcomprises determining reference coordinates for points A, B, and C ofthe vessel, generating a signal representative of velocity at a point onthe vessel, generating a heave rate signal based upon the velocitysignal, and generating a heave signal based upon the heave rate signal.

[0005] In another aspect, an inertial navigation system (INS) fordetermining a heave and a heave rate for a vessel is provided. The INSsystem comprises, in an exemplary embodiment, a main unit, a userinterface, a global positioning satellite (GPS) receiver, and a sensorunit. The sensor unit comprises an inertial sensor assembly located at apoint A of the vessel. The main unit comprises interfaces forcommunication of inertial navigation information to other controlsystems on the vessel. The INS is configured to determine referencecoordinates for points A, B, and C of the vessel, a heave, and a heaverate for the vessel.

[0006] In a further aspect, a filter is provided which comprises a firststage and a second stage. The first stage is configured to generate anoutput signal representative of a heave rate based, at least in part ona velocity signal input. The second stage is configured to generate anoutput signal representative of heave by filtering the heave rate signaloutput by the first stage.

BRIEF DESCRIPTION OF THE DRAWINGS

[0007]FIG. 1 is a block diagram of an inertial navigation system.

[0008]FIG. 2 is a diagram illustrating factors utilized in determiningheave and heave rate.

[0009]FIG. 3 is a block diagram of a filter.

[0010]FIG. 4 is a block diagram of a filter utilizing a different inputthan the filter of FIG. 3.

[0011]FIG. 5 is a detailed block diagram of two filters, illustratingtransfer functions for each filter.

DETAILED DESCRIPTION OF THE INVENTION

[0012] The features and principles are described herein relative to anexemplary embodiment thereof. It will be apparent to those skilled inthe art that numerous variations and modifications may be made to theexemplary embodiment without departing from the spirit and scope. Thesystems and methods are not limited to the specific embodimentsdescribed herein. Components of each system and method can be practicedindependent and separate from other components and methods. Each systemand method also can be used in combination with other components,systems, and methods.

[0013] A system that provides an accurate representation of heave andheave rate of a vessel, as measured by a strapdown inertial system(INS), is described in detail below. Generally, the system utilizesdigital band pass filtering to remove long-term bias errors due tospecific force measurement and modeling of local gravity. In addition,the system performs lever-arm corrections to define velocity at pointsof the vessel other than points at which sensors of the INS are located,using measured values of angular velocity of the vessel. As used herein,the term vessel means any structure which may incorporate a system fordetermining heave and heave rate including, but not limited to, ships,boats, submarines, submersibles, marine vessels, unmanned watervehicles, military vessels, amphibious aircraft, torpedoes, and anyother watercraft in a marine environment.

[0014]FIG. 1 is a block diagram of an inertial navigation system (INS)10 for marine use. System 10 includes a main unit 12 to which isconnected a power supply 14. Also connected to main unit 12 are a userinterface unit 16, a global positioning satellite (GPS) receiver 18, anda sensor unit 20 which includes an inertial sensor assembly (ISA) 22.Main unit 12 is configured with interfaces 24 to provide inertialnavigation information to, for example, a main navigation system withina vessel, or other systems within a vessel that may utilize navigationinformation. ISA 22, within sensor unit 20 incorporates gyroscopes andaccelerometers (not shown) and sensor unit 20 includes gyroscopeelectronics. In one embodiment, the gyroscope is a laser gyro. Userinterface unit 16 provides a user interface, including a display andkeypad (neither shown) for display of inertial information andalteration of display modes and operating parameters.

[0015]FIG. 2 is a diagram of a vessel 40 illustrating a location of ISA22, herein denoted as point A. Point B is referred to as a heave zeroreference point. In one embodiment, point B is chosen to be at a meanwater line and at a mean center of rotation of vessel 40. One purposefor designating point B is to compensate for a static tilt of vessel 40.Another purpose for the designation is compensation for how high or lowvessel 40 is riding in water 42 due to buoyancy. Point C is a point atwhich it is desirable to know heave and heave rate for vessel 40, forexample, where a winching system is connected to a submersible (notshown).

[0016] To determine heave and heave rate for vessel 40, referencecoordinates are defined. Specifically, body frame coordinates aredefined with respect to a bow, a starboard side and a belly of vessel40. Local level frame coordinates are defined as north, east, and down.C_(B) ^(L) is a transformation direction cosine matrix that redefines avector expressed in body frame coordinates into a vector expressed inlocal level coordinates. For transformation matrices, the subscriptindicates the frame the vector is to be taken from, and the superscriptindicates the frame the vector is to be taken to. Therefore referring tothe above transformation direction cosine matrix, the “B” refers to bodyframe coordinate vectors, and “L” refers to local level frame vectors.

[0017] A lever arm from point A to point C is given as R_(AC), and alever arm from point A to point B is given as R_(AB). Lever arms, asknown in the art, are vectors representing distance and direction. Leverarms herein are defined in body frame coordinates. Referring to FIG. 2,a lever arm from point B to point C can be described as:

R _(BC) =R _(AC) −R _(AB)  Eq (1)

[0018] Velocity of vessel 40 is known at point A, based upon data fromINS 10. Therefore velocity at point C is calculated as:

V _(C) ^(L) =V _(A) ^(L) +C _(B) ^(L)(ω^(B) ×R _(AC))  Eq (2)

[0019] Superscripts on the above variables indicate which frame, bodycoordinates or local level coordinates, the variable is defined in.Therefore, V_(A) ^(L) is the local level velocity at point A, V_(C) ^(L)is the local level velocity at point C, and ω^(B) is the body frameangular rate of vessel 40.

[0020] In one embodiment, inertial navigation system (INS) 10 includes afilter module that provides heave and heave rate by filtering analtitude channel of INS 10. In at least some known applications, thealtitude channel of INS 10 is augmented with external measurements froma global positioning satellite (GPS), for example. Other applicationsprovide augmentations for the altitude channel of INS 10 from one ormore of GPS, pressure altitude, or simply loosely slaving the altitudechannel of INS 10 to zero (i.e. sea-level altitude).

[0021] To remove bias errors, such as bias errors from altitude andaltitude rate, high pass filtering is used to arrive at signalsrepresentative of heave and heave rate, thereby avoiding the errorsassociated with GPS, baro-altimeter or other measurement systems.

[0022] Referring to FIG. 3 a filter module within INS 10 is shown. Athird element of V_(C) ^(L), referred to herein as V_(C) ^(L)(3), isinput into a high pass filter 60 to remove any DC or very low frequencycomponents. Providing an input 62 of V_(C) ^(L)(3) to filter 60 resultsin an output 64 of {circumflex over (V)}_(C) ^(L)(3). Output 64, thesignal that results from filtering through filter 60, is the heave rateof vessel 40. Heave rate is therefore stated mathematically as follows:

HeaveRate={circumflex over (V)}_(C) ^(L)(3)  Eq (3)

[0023] In one embodiment of system 10 (shown in FIG. 1), a 200 Hzattitude rate is utilized in conjunction with Eq (2) which contains aninherent white noise variance of 6.0e-6 (radians)²/(second)².

[0024] A change in vertical position of vessel 40, that is, heave orheave distance, is then calculated by integration of the heave ratesignal. An example calculation for heave signal, utilizing a heave ratesignal follows: $\begin{matrix}{{Heave} = {{\sum{{{\hat{V}}_{C}^{L}(3)}\Delta \quad t}} + {C_{B}^{L}{{R_{BC}(3)}_{t_{0}}.}}}} & {{Eq}\quad (4)}\end{matrix}$

[0025] In Eq (4), the first term, ∑V̂_(C)^(L)(3)Δ  t,

[0026] represents heave rate integration, and the second term, C_(B)^(L)R_(BC)(3)_(t) ₀ , is the initial condition on the integration. Theinitial condition represents a local level vertical distance betweenpoint B and point C, which was valid at time to, a starting time of theintegration.

[0027] Since white noise is inherent in the attitude rate, as determinedby INS 10, output 64 of high pass filter 60 may be unpredictable wheninput 62 contains such white noise. It is therefore desirable that thedata to be filtered in a for determination of heave and heave rate be asaccurate as possible. In the embodiments described herein, the mostaccurate data is the data for point A. Therefore, an alternativeembodiment of a heave algorithm implemented within a filter module inINS 10 is provided which results in a calculation for heave and heaverate with results similar to those described above. This alternativeembodiment, illustrated in FIG. 4, introduces less noise.

[0028] Referring specifically to FIG. 4, high pass filter 70 receives alocal level vertical velocity signal at an input 72, the local levelvertical velocity is provided from INS 10 at point A, as represented bythe input variable in the Figure. In one specific embodiment, locallevel vertical velocity is a signal provided at 200 Hz. In such anembodiment, a heave rate signal is calculated as: $\begin{matrix}{{V_{C}^{L}(3)} = {{{\hat{V}}_{A}^{L}(3)} + {\left\lbrack {C_{B}^{L}\left( {\omega^{B} \times R_{A\quad C}} \right)} \right\rbrack {(3)_{t_{k}}.}}}} & {{Eq}\quad (5)}\end{matrix}$

[0029] The subscript t_(k) in the second term of Eq (5) indicates thatthe term is updated every iteration. Utilizing Eq (4), heave is thencalculated as an integration of Eq (5), which follows: $\begin{matrix}{{\sum{{V_{C}^{L}(3)}\Delta \quad t}} = {{\sum{{{\hat{V}}_{A}^{L}(3)}\Delta \quad t}} + \left\lbrack {C_{B}^{L}{R_{BA}(3)}_{t_{k}}} \right\rbrack + {\left\lbrack {C_{B}^{L}R_{A\quad C}} \right\rbrack (3)_{t_{k}}}}} & {{Eq}\quad (6)}\end{matrix}$

[0030] The subscript t k in the second term and in the third termindicates that the term is updated every iteration. Equations 5 and 6provide signals for heave and heave rate which are similar to thoseprovided by filter 60 (shown in FIG. 3), but with less potential fornoise and unpredictable behavior in the filtering process.

[0031]FIG. 5 illustrates a specific embodiment of a filter 100implemented within INS 10. Filter 100 determines heave rate utilizingvertical velocity, and a filter 102 which determines a filtered heavedistance utilizing heave rate. Referring specifically to filter 100,vertical velocity is received at input 104 and heave rate is placed atoutput 106. Referring specifically to elements within filter 100, afirst filter element 108, which has a transfer function of(1+z⁻¹)(1+z⁻²)/4, receives vertical velocity as input. An output offirst filter element 108 is input to both of a positive input of asubtraction element 110 and an input of a second delay element 112.Second delay element 112 has a transfer function of(1+z⁻⁴)(1+z⁻⁸)/(15,500−30,750z⁻¹⁶+15,254z⁻³²). An output of element 112is input to a third delay element 114. Third delay element 114 has atransfer function of (1+z⁻¹⁶)/(125−123z⁻¹⁶). An output of third delayelement 114 is input to a negative input of subtraction element 110. Anoutput of subtraction element 110 is provided to an input of a fourthdelay element 116 that has a transfer function of 4/(33−48z⁻¹+19z⁻²).Output of fourth delay element 116 is output 106, and therefore theheave rate of vessel 40.

[0032] Output of fourth delay element 116, or heave rate is input tofilter 102, at a first delay element 118. First delay element 118 has atransfer function of ({fraction (1/200)})/(1−z⁻¹). Output of first delayelement 118 is heave 120. Heave 120 is input to a positive input ofsubtraction element 122 and an input of a second delay element 124.Second delay element 124 has a transfer function of(1+z⁻⁴)(1+z⁻⁸)/(15,500−30,750z⁻¹⁶+15,254z⁻³²). An output of second delayelement 124 is input to a third delay element 126. Third delay element126 has a transfer function of (1+z⁻¹⁶)/(125−123z⁻¹⁶). An output ofthird delay element 126 is input to a negative input of subtractionelement 122. An output of subtraction element 122 is the filtered heavedistance of vessel 40.

[0033] Utilization of filters 100 and 102 within INS 10, in oneembodiment, provide signals with an accuracy of at least 0.1 meter forheave and 0.1 meter per second for heave rate. Filters 100 and 102further provide a pass band of 0.1 to 2 Hertz and a group delay of lessthan 25 milliseconds. In the embodiments described, heave has a range ofplus to minus four meters.

[0034] A particular approach is provided in the elements of filters 100and 102. An upper edge, of the filters, out to 2 Hertz and above isconfigured with unity gain, within 0.015 dB, and with anti-aliasingzeros for either 50 or 200 sample per second outputs. Instead of twozeros at 100 Hertz, as is typical in certain known bi-lineartransformation second order complex-pole filters, one of the zeros ismoved to 50 Hertz. Filters 100 and 102 therefore do not exhibit aliasingerrors, due to near 50 and 100 Hertz content, when filter output isutilized at 50 samples per second.

[0035] A low-frequency end of filters 100 and 102 is configured as ahigh pass filter and has unity gain, within about 0.20 dB, at 0.1 Hertz.A ratio of filter break frequency to sampling frequency, raised to thepower equal to the number of filter poles, causes small differences infilter coefficients which thereby locates the poles of the filters.Although data from filters 100 and 102 is available at a high samplingrate, the poles are realized in a computation repeated at a slow rate.The purpose of the low-frequency end of filters 100 and 102 is toextract an average value for the signals.

[0036] The reduction in processing rate is therefore done in a mannerthat attenuates signal content that can alias to a near-zero frequency,before the reduction in frequency occurs. An analysis of such a processshows that content near the zeros is removed before an aliasing to anear-zero frequency occurs with the reduction in sampling rate. Furtherreduction in effect of small differences of large coefficients is toimplement the overall transfer function as a cascade of first andsecond-order sections, rather than implementing as higher ordersections.

[0037] While the invention has been described in terms of variousspecific embodiments, those skilled in the art will recognize that theinvention can be practiced with modification within the spirit and scopeof the claims.

What is claimed is:
 1. A method for determining heave and heave rate ata point C on a vessel, the vessel having an inertial navigation system(INS), sensors for the INS being located at a point A on the vessel, thevessel having a zero heave reference point B, said method comprising:determining reference coordinates for points A, B, and C; generating avelocity signal representative of a velocity at a point on the vessel;generating a heave rate signal based on the velocity signal; andgenerating a heave signal based upon the heave rate signal.
 2. A methodaccording to claim 1 wherein determining reference coordinatescomprises: determining body coordinates for the vessel at points A, B,and C; and determining local level coordinates for the vessel at pointsA, B, and C.
 3. A method according to claim 2 wherein, for point C ofthe vessel, determining a velocity of the vessel comprises: defining atransformation direction cosine matrix C_(B) ^(L), which redefines avector expressed in body frame coordinates into a vector expressed inlocal level coordinates; determining a lever arm from point A to pointC, R_(AC); and calculating a velocity at point C, using a known velocityat point A according toV_(C)^(L) = V_(A)^(L) + C_(B)^(L)(ω^(B) × R_(A  C)), where  V_(A)^(L)

is the local level velocity at point A, V_(C) ^(L) is the local levelvelocity at point C, and ω^(B) is the body frame angular rate of thevessel.
 4. A method according to claim 3 wherein filtering the velocitysignal to provide a heave rate signal comprises filtering a thirdelement of V_(C) ^(L), V_(C) ^(L)(3), with a high pass filter todetermine a heave rate, {circumflex over (V)}_(C) ^(L)(3).
 5. A methodaccording to claim 4 wherein integrating the heave rate signal todetermine a heave comprises calculating a heave according toHeave = ∑V̂_(C)^(L)(3)Δ  t + C_(B)^(L)R_(BC)(3)_(t₀), where  V̂_(C)^(L)(3)Δ  t,

represents heave rate integration, and C_(B) ^(L)R_(BC)(3)_(t) ₀ , is aninitial condition on the integration, the initial condition being alocal level vertical distance between point B and point C, which wasvalid at time t₀, the starting time of the integration.
 6. A methodaccording to claim 2 further comprising: defining a transformationdirection cosine matrix C_(B) ^(L), which redefines a vector expressedin body frame coordinates into a vector expressed in local levelcoordinates; and determining a lever arm from point A to point C,R_(AC).
 7. A method according to claim 6, wherein velocity of the vesselis known at point A, and wherein filtering the velocity signal toprovide a heave rate signal comprises determining heave rate of thevessel according to V_(C) ^(L)(3)={circumflex over (V)}_(A)^(L)(3)+[C_(B) ^(L)(ω^(B)×R_(AC))](3)_(t) _(k) , where {circumflex over(V)}_(A) ^(L)(3) is a filtered third element of local level velocity atpoint A, V_(C) ^(L)(3) is the third element of local level velocity atpoint C, ω^(B) is the body frame angular rate of the vessel, and thesubscript t_(k) indicates that the term is updated every iteration.
 8. Amethod according to claim 7 wherein integrating the heave rate signal todetermine a heave comprises determining a heave according to∑V_(C)^(L)(3)Δ  t = ∑V̂_(A)^(L)(3)Δ  t + [C_(B)^(L)R_(BA)(3)_(t_(k))] + [C_(B)^(L)R_(A  C)](3)_(t_(k)),

where the subscript t_(k) indicates that the term is updated at everyiteration.
 9. A method according to claim 1 wherein filtering thevelocity signal to provide a heave rate signal comprises filtering thevelocity signal according to the transfer function[(1+z⁻¹)(1+z⁻²)/4]×[4/(33−48z⁻¹+19z⁻²)]×[1−((1+z⁻⁴)(1+z⁻⁸)(1+z⁻¹⁶)/(15,500−30,750z⁻¹⁶+15,254z⁻³²)(125−123z⁻¹⁶))].
 10. A methodaccording to claim 1 wherein integrating the heave rate signal toprovide a heave comprises filtering the heave rate signal according tothe transfer function ({fraction (1/200)})/(1−z⁻¹).
 11. A methodaccording to claim 10 further comprising filtering the heave accordingto 1−[(1+z⁻⁴)(1+z⁻⁸)(1+z⁻¹⁶)/(15,500−30,750z⁻¹⁶+15,254z⁻³²)(125−123z⁻¹⁶)].
 12. An inertial navigation system (INS) fordetermining a heave and a heave rate for a vessel, said systemcomprising: a main unit; a user interface unit; a global positioningsatellite (GPS) receiver; and a sensor unit comprising an inertialsensor assembly (ISA), said ISA located at a point A of the vessel, thevessel having a zero heave reference point B, said main unit comprisinginterfaces for communication of inertial navigation information to othersystem on the vessel, said user interface unit comprising a display andkeypad, said system configured with reference coordinates for points A,B, and C, point C being a location at which it is desired to be providedheave and heave rate data, said system configured to determine avelocity of the vessel at a point on the vessel, velocity beingrepresented by a signal, filter the velocity signal to provide a heaverate signal, and integrate the heave rate signal to determine a heave.13. An INS system according to claim 12 wherein the referencecoordinates for points A, B, and C of the vessel comprise: bodycoordinates for the vessel at points A, B, and C; and local levelcoordinates for the vessel at points A, B, and C.
 14. An INS systemaccording to claim 13 wherein to determine a velocity for point C of thevessel said system is configured to: define a transformation directioncosine matrix C_(B) ^(L), which redefines a vector expressed in bodyframe coordinates into a vector expressed in local level coordinates;determine a lever arm from point A to point C, R_(AC); and calculate avelocity signal at point C, using a known velocity at point A accordingto V_(C) ^(L)=V^(L) _(A)+C_(B) ^(L)(ω^(B)×R_(AC)), where V_(A) ^(L) isthe local level velocity at point A, V_(C) ^(L) is the local levelvelocity at point C, and ω^(B) is the body frame angular rate of thevessel.
 15. An INS system according to claim 14 wherein to filter thevelocity signal to provide a heave rate signal, said system isconfigured to filter a third element of V_(C) ^(L), V_(C) ^(L)(3), witha high pass filter to determine a heave rate signal, {circumflex over(V)}_(C) ^(L)(3).
 16. An INS system according to claim 15 wherein tointegrate the heave rate signal to determine a heave, said system isconfigured to calculate a heave according toHeave = ∑V̂_(C)^(L)(3)Δ  t + C_(B)^(L)R_(BC)(3)_(t₀), where  V̂_(C)^(L)(3)Δ  t,

represents heave rate integration, and C_(B) ^(L)R_(BC)(3)_(t) ₀ , is aninitial condition on the integration, the initial condition being alocal level vertical distance between point B and point C, which wasvalid at time t₀, the starting time of the integration.
 17. An INSsystem according to claim 13 wherein said system is configured to:define a transformation direction cosine matrix C_(B) ^(L), whichredefines a vector expressed in body frame coordinates into a vectorexpressed in local level coordinates; and determine a lever arm frompoint A to point C, R_(AC).
 18. An INS system according to claim 17,wherein said system is provided a velocity of the vessel at point A, andwherein to filter the velocity signal to provide a heave rate signal,said system is configured to determine heave rate of the vesselaccording to V_(C) ^(L)(3)={circumflex over (V)}_(A) ^(L)(3)+[C_(B)^(L)(ω^(B)×R_(AC))](3)_(t) _(k) , where {circumflex over (V)}_(A)^(L)(3) is a filtered third element of local level velocity at point A,V_(C) ^(L)(3) is the third element of local level velocity at point C,ω^(B) is the body frame angular rate of the vessel, and the subscriptt_(k) indicates that the term is updated every iteration.
 19. An INSsystem according to claim 18 wherein to integrate the heave rate signalto determine a heave, said system is configured to determine a heaveaccording to∑V_(C)^(L)(3)Δ  t = ∑V̂_(A)^(L)(3)Δ  t + [C_(B)^(L)R_(BA)(3)_(t₀)] + [C_(B)^(L)R_(A  C)](3)_(t_(k)),

where the subscript t₀ indicates that the term is an initial conditionand is only calculated once at the start of integration.
 20. A filtercomprising: a first stage configured to provide a heave rate signalbased upon a velocity signal input; and a second stage configured toprovide a heave signal by filtering the heave rate signal output by saidfirst stage.
 21. A filter according to claim 20 wherein said first stageis configured to filter the velocity signal input according to thetransfer function[(1+z⁻¹)(1+z⁻²)/4]×[4/(33−48z⁻¹+19z⁻²)]×[1−((1+z⁻⁴)(1+z⁻⁸)(1+z⁻¹⁶)/(15,500−30,750z⁻¹⁶+15,254z⁻³²)(125−123z⁻¹⁶))],and said second stage is configured to filter the heave rate signalaccording to the transfer function ({fraction (1/200)})/(1−z⁻¹).
 22. Afilter according to claim 21 wherein said second stage is configured tofilter the heave signal according to according to the transfer function1−[(1+z⁻⁴)(1+z⁻⁸)(1+z⁻¹⁶)/(15,500−30,750z⁻¹⁶+15,254z⁻³²)(125−123Z⁻¹⁶)],providing a filtered heave distance signal.